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ABSTRACT 



Although it has important ramifications for both the formation of star clusters and their subsequent dynamical evolution, rotation 
remains a largely unexplored characteristic of young star clusters (few Myr). Using multi-epoch spectroscopic data of the inner 
regions of 30 Doradus in the Large Magellanic Cloud (LMC) obtained as part of the VLT-FLAMES Tarantula Survey, we search 
for rotation of the young massive cluster R136. From the radial velocities of 36 apparently single O-type stars within a projected 
radius of 10 pc from the centre of the cluster, we find evidence, at the 95% confidence level, for rotation of the cluster as a whole. We 
use a maximum likelihood method to fit simple rotation curves to our data and find a typical rotational velocity of ~3 km s _1 . When 
compared to the low velocity dispersion of R136, our result suggests that star clusters may form with at least ~20% of the kinetic 
energy in rotation. 

Key words, galaxies: star clusters: individual (R136) - Magellanic Clouds - stars: kinematics and dynamics - globular clusters 
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1. Introduction 

Despite their spherical shape, Milky Way globular clusters 
(GCs) rotate with amplitudes up to half th e ID velocity dis- 
persio n (0 < V mt sin i /o"id < 0.5, e.g. iMevlan & HeggieL 
[l997h . so their amount of rotational energy is typically not dom- 
inant but also not negligible. Numerical studies have shown that 
this rotation has an important influence on star clusters by ac- 
celerating their dynamical evolution, for example by speeding 
up the collapse of the core through the gravogyro instability 
or by significantly i ncreasing the escape rate for clusters in a 
strong tidal field ( e.g. lEinsel & Spurzemlll999tlKim et al.Ll2002t 
lEmstetaUl2007n . 

Most of the rotation signatures are found through radial 
velocity (RV) studies, but rotation has also been confirmed 
in the plane of th e sky for to Centauri and 47 Tucanae 



sky 

(Ivan Leeuwen et al.ll2000r . lA"nderson & van der Mare]Ll2010l re- 
spectively). RV studies are now able to measure rotational am- 
plitudes in GCs below 1 kms 1 and despite these precise mea- 
surements rotation has not been detected in some clusters (e.g. 



* Based on observations collected at the European Southern 
Observatory under program ID 182.D-0222 



Lane et al although this could also be an inclination ef- 

fect. 

It is unclear what the origin of the rotation is in some of these 
old clusters. It could be the result of the merging of two clusters 
(Bau mgardt et al., 2003) or imprinted during the formation pro- 
cess. Observational input is now getting sufficiently abundant 
to look for correlati ons between ro t ationa l amplitude and other 
cluster properties. Bellazzin i et alj d2012l) report a correlation 
between horizontal branch (HB) morphology and V 10t sin / /<x ID 
in a sample of 20 globular clusters. Given that metallicity is the 
first parameter determining HB morphology, this in turn sug- 
gests a correlation between V rot sin« /<Xi D and metallicity, such 
that clusters with higher metallicity have greater fractions of en- 
ergy in rotation. Since a higher metallicity in a gas implies a 
higher efficiency of energy dissipation by atomic transitions, this 
then hints at a signifi cant ro l e of dissipation in the process of 
cluster formation (e.g. Bekki l201Q|) . Rotation may therefore be 
intimately linked to the formation of clusters. 

Little is known about rotation in young clusters, par- 
tially because it is very challenging to measure accurately 
am given the high multipl i city fraction of mas sive stars (e.g. 
Henault- Brunet et al.Ll2T3l2l) . lFrenk & Falll(11982l) found an age- 
ellipticity relation for clusters in the Large Magellanic Cloud 
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(LMC) with older clusters presenting less elongated shapes, 
which was interpreted as internal evolution erasing any asym- 
metry stemming fro m the violent rel axation process of the for- 
mation (although see Goodwill 1 19971) . However, rotation and el- 
lipticity are not nece ssarily equivalent. Ellipticity can be due to 
rotation (e.g. u > Cen; | Mevlan & Mayorlll986 ) but also to veloc- 
ity anisotropy (Step hens et all l2006tlHenonL| 19731), and ro tating 
clusters can be spherical ( Lvnden-Belll [i960; Mez a, 120021) . 

Marginal evidence for rotation was found for the young 
(few 100 Myrs) Galactic clust er GLIMPSE-CO 1 w ith an ampli- 
tude of V rot sin/ /cr 1D ^ 0.2 (Dav ies et al rotational 
signal in t he RVs was also de tected in the ~100 Myr cluster 
NGC 1 866 dFischer et all 1 19921) an d in the -50 Myr binary clus- 
ter NGC 1850 dFischer et all 1 19931) . both in the LMC. To really 
confirm whether clusters form with a significant amount of an- 
gular momentum, we need to look for an even younger clus- 
ter. The young massive cluster (YMC) R136 in the 30 Doradus 
star forming region in the LMC is an ideal ta rget to establish 
this. W ith an estimated mass of about 10 5 M o (lAndersen et all 
2009) and its sub solar metallicity, it may at some stage resem- 
ble a typical metal-rich GC as we find th em in the Milky Way 
Bulge. With an age of l ess than 2 Myr (Ide Koter et all 1 19981: 
iMassev & Hunter! 119981 ICrowther etafl 120101) , it is "so young 
that any rotation needs to be attributed to the formation process, 
be it from merging of sub-clusters or directly from the angu- 
lar momentum of the progenitor cloud. A rough estimate of the 
half-mass relaxation time (f r h) of R 136 can be obtained by as- 
suming N = 10 5 stars and a half-mass radius of 2.27 pc, which 
is found from multiplying the half-li ght rad ius of 1.7pc (e.g. 
IHenault-Brunetet all 1201 2|) by 4/3 (ISpitzerl 119871) . Following 
the formula of Snit zer & Hard (Il97ll) . we obtain f,h - 366 Myr, 
so relaxation would not have had time to erase the original sig- 
nature of rotation. 

Here we report on evidence for rotation of R136 deduced 
from RV measurements of massive stars obtained as part o f the 
VLT-FLAMES Tarantula Survey (VFTS: lEvans^tini20Tlh . We 
briefly present the data in Sect.|2]and describe our analysis of the 
rotational signature in Sect. [3] We discuss the implications of the 
rotation of R 136 for cluster evolution in Sect. [4] and present our 
conclusions in Sect. [5] 
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Fig. 1. Illustration of the positions and RVs of the stars consid- 
ered in this study. Symbol sizes denote the magnitude of the 
stellar velocities with respect to the average cluster velocity. The 
solid, dotted, and dashed lines correspond to the optimal rotation 
axis determined for models with a constant rotational velocity, a 
constant rotation rate, and a more realistic rotation curve (see 
Sect. [3]), respectively. 



all between 5 and lOpc except two between 4 and 5pc. The 
identification numbers, coordinate s, spectral types, and RVs o f 
these stars are listed in Table 3 of Henault-Brun et et al ] (120121) . 
In Fig.Q] we schematically present their positions and distribu- 
tion of RVs with respect to the mean RV of the cluster. 



3. Analysis 

To look for the signature of rotation, we fit rotating models to 
our set of measure d RVs by maximi zing the logarithm of the 
likelihood function (Bevingtorl ll969l) 



2. Data 

The data used in this work consist of RV measurements and 
their uncertainties for 36 apparently single O-type stars within 
a projected radius of lOpc from the centre of R136 (adopted 
here as the position of the star R136-al: a = 5 h 38 m 42 s . 39, 6 = 
-69°06'02."91, J2000). These are based on observations from 
the VFTS, for which at least five epochs were obtained for five 
different pointings of the FLAMES-ARGUS integral-field unit 
in the central arcminute of 30 Dor, in addition to nine config- 
urations of the Medusa fibre-feed to the Giraffe spectrograph 
in the surround ing Tarantu la Nebula over a 25 arcmin diameter 
field-of-view (Evans et ail 1201 ll) . The RV and variability analy- 
sis, based on Gaussian fitting of selected helium stella r absorp- 
tion lines, is presen ted in iHenault-Brunet et all (1201 2l) for the 
ARGUS data and in lSana et all d2012l) for the Medusa data. We 
retain here only the RVs of stars showing no significant variabil- 
ity throughout all epochs, and we apply a 3 <x clipping centered 
on the mean velocity of the cluster, with <x = 6kms _1 , the ob- 
served line-of-sight velocity dispersion of the cluster (i.e. before 
correction of the velocity di spersion for undetected binaries, see 
IHenault-Brunetet all |2012|) . This yields a total of 16 ARGUS 
sources, all within 5 pc from the centre, and 20 Medusa sources, 



M = \n£ = In 



(i) 



Given our relatively small data set, we only consider three simple 
models: constant rotational velocity amplitude, constant rotation 
rate, and finally a more realistic model with solid body rotation 
in the inner parts peaking near 2 half-mass radii followed by a 
decline. For these models, the probability density function Pi of 
a measurement RV,- with uncertainty cr,- given a position angle 
PAo for the rotation axis and a line-of-sight velocity dispersion 
cr 1D can be written as 



O 2 + (T 2 lD 



: exp 



RV, - V rot sin i 



2 
ID 



2 



(2) 



where V rot sin i is a constant for a model with fixed rotational 
velocity. For a model with constant rotation rate, V rot sin i de- 
pends on the rotation rate (£2) and the distance from the rota- 
tion axis (Xj) such that [V rot sin/]j = Q Xj. For the physically 
motivated model, we adopt a function of the form [V 10t sin i]j = 
A/2-Xj/(l +(Xj/4) 2 ), where Xj is the distance to the rotation axis 
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in pc and A is the maximum rotational velocity which we assume 
is at 4pc for simplicity (the half-mass radius of R 136 is about 
2pc). This rotation curve captures the general behaviou r seen 
in simulations of rotating clusters (e.g. iKim et all 12002). Note 
that b ecause the velocity d i spersi on profile of R 136 is relatively 
flat dHenault-Brunet et al.L 120121) . this function also describes 
V r ot/cr. We define the position angle (with respect to the centre of 
the cluster) as increasing anti-clockwise in the plane of the sky 
from North (PA= 0°) towards East (PA= 90°). We adopt neg- 
ative rotational velocities for position angles between PAo and 
PAo+180°. For simplicity, we assume that cr 1D is constant across 
the radius range considered. This cr 1D is largely due to cluster 
dynamics, but contains a small contributi on from the orbital mo- 
tions of undetected binaries (1-2 km s" 1 : IHenault-Brunet et ail 
1201 2|) . effectively adding noise to the rotation signature we are 
trying to detect. 

The values of the parameters that maximize the likeli- 
hood function are PA =57±15°, V rot sin/=2.9±1.0kms _1 and 
crirj=5.9±0.8kms~' for the constant rotational velocity model, 
PA =44+18°, Q=0.75+0.22Myr 1 and cr 1D =6.0+l.l kms" 1 
for the constant rotation rate model, and PAo=36+14°, 
A=3.7+1.3kms~ 1 and cr 1D =5.9+0.8kms _1 for the final model. 
The best-fit rotation curves are shown in Fig. [2] The uncertain- 
ties were estimated using Monte Carlo calculations on simulated 
data sets comparable to our measured data (i.e. with the same 
rotational signal). We also obtained consistent uncertainties by 
considering the parameter change neces sary to decrease M by 
AM = | from its value at the maximum dBevingtonl[T969l) . The 
maximum value of M is slightly higher for the constant rota- 
tional velocity model compared to the two other models, but 
likelihood ratio tests showed that the difference is not significant. 
Thus, we currently cannot favour one model over the others. 

To establish the significance of the detected rotational signal, 
we performed Monte Carlo simulations and applied the maxi- 
mum likelihood method above to 10000 random distributions 
of velocities (i.e. non-rotating systems). The adopted size and 
spatial distribution of the simulated populations were taken to 
be the same as the observed sample in order to be sensitive to 
possible biases introduced by our non-uniform spatial sampling. 
The velocities were drawn from a Gaussian distribution with <x 
determined by the observed line-of-sight velocity dispersion. For 
each simulated star, observational noise was added based on the 
RV uncertainty of the observed star at the same location, again to 
take into account the biases introduced by our non-uniform data 
set. The distribution of Q., V rot sin i and A from 10 000 such sim- 
ulations is shown in Fig. [2] The distributions do not peak at zero 
because the limited number of stars and the measurement noise 
generally result in a non-zero amplitude. Our best-fit values are 
located in the right tail of these distributions and the correspond- 
ing confidence level of the detection is around 95% for all three 
models. 

The analysis outlined above was also performed on a sub- 
sample from which supergiant candidates and stars with possible 
composite spectra were excluded as their RVs might be inaccu- 
rate (see Hen ault-Brunet et a P. 120121) . Very similar results were 
obtained, with best-fit parameters all within 4% of those previ- 
ously determined and a confidence level about 1 % lower because 
of the smaller number of measurements. As an additional check, 
we ran the analysis on different subsamples of appare ntly single 
O-typ e stars from VFTS using the measurements of ISana et alj 
(1201 2l) . including stars farther away from R136 than in our main 
sample. The analysis of these subsamples suggests that the rota- 
tion signal is dominated by the stars in the inner regions and does 
not result from a velocity gradient across the field on a larger 



scale in the surrounding OB association. For example, when con- 
sidering the sample of stars between 10 and 20 pc from the cen- 
tre (37 stars), we find that V rot sin/ goes down to 1.8kms~' for 
a constant rotational velocity model and the confidence level of 
the rotational signature is only 51%. Note that stars that are part 
of the surrounding 30 Doradus region and not members of R 136 
could contaminate our sample. If we assum e that the outer co m- 
ponent of the double-co mponent EFF fit (Elson et alj 119871) to 
the light profile of R 136 (Ma ckev & Gilmorell2003h is due to the 
surrounding association, then we may expect these stars to con- 
tribute to > 50% of the sample beyond 5 pc from the centre and 
to dilute the rotational signal in the outer regions of the cluster. 

We also looked for the signature of rotation with a method 
commonly used in studies of globular clusters. This method con- 
sists of dividing the sample with a line passing through the centre 
of the cluster and computing the difference in the mean RV be- 
tween the two subsamples of stars as a function of the position 
angle of the dividing line (e.g. lCote et al.L[l995l : lBellazzini et all 
120121) . Results very similar to those reported above were ob- 
tained for the position of the rotation axis and the rotational am- 
plitude, but the maximum likelihood method has the advantage 
of directly comparing the data to models without binning. 



4. Discussion 

Given that <x 1D 5 ± 1 kms" 1 for R136 after correcting 
for undetected binaries (IHenault-Brunet et all 1201 2l) and that 
the mean rotational velocity in the radius range considered is 
=; 3 ± 1 km s _1 for the three best-fit rotation curves shown above, 
our analysis implies that V rot sin/ /cr 1D 0.6 + 0.3, which is 
somewhat larger than what is typic ally observed in globular clus- 
ters (e.g. lMevlan & HeggieL[l9 97). We have to keep in mind that 
y ro t sin / is itself a projected quantity, and so is cr 1D . We would 
like to know how these quantities relate to the 3D quantities (v^) 
and v rms . Here (v#) is the mean tangential velocity component in 
the (p direction, used as the rotational component of the velocity 
vector in spherical coordinates for a cluster rotating about the z 
axis, and v rms is the root-mean square of the velocities. 

We can estimate that V rot /cr\D is typically about 10% larger 
than (v^/Vrms. This can be understood as follows. For a circu- 
lar orbit with unit velocity and the line of sight along the orbital 
plane (sin i = 1) the mean RV component is 2/n. For an isotropic 
velocity dispersion the observed component is <x 1D = v rms / y3. 
This means that [^ ro t/o"iD]/[<v^>/v rms ] - 2 V3/7T 1.1. If 
sin/ = 1 then we would have <v^)/v rms 0.9V rot /cr 1D and a 
lower limit on (v^,)/v rms of 0.5±0.2. A ratio of (v0)/v rms = 0.5 
implies a ratio of kinetic energy in rotation over kinetic energy 
in random motions of 0.25, which in turn implies that at least 
20% of the total kinetic energy is in rotation. If instead we as- 
sume that sin/ = 2/n (the average of sin / assuming a random 
distribution of inclination angles), we have (v^)/v ims = 0.9 and 
45% of the total kinetic energy in rotation. The approximate cri- 
terion for stability against non-axisy mmetric perturbations pro- 
vided by lOstriker & Peebles! d 1973b states that the ratio of the 
rotational kinetic energy over the potential energy should not 
exceed \T mt /W\ < 0.14. If we assume virial equilibrium, then 
T = -0.5W, so for T TOt = 0.20 T and T Iot = 0.45 T this implies 
|7rot/Wl = 0.1 and 0.225, respectively. Although not a rigorous 
test, this suggests that a low inclination is perhaps more reason- 
able. 

If the cluster was flattened by rotation (although recall our 
previous words of caution about ellipticity and rotation), we 
would expect to see peaks in the azimuthal density profile, with 
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Fig. 2. Top: RV in the system of the cluster as a function of distance X from the optimal rotation axis for models with constant 
rotation rate (left), constant rotational velocity (centre), and the more realistic rotation curve discussed above (right). The best-fit 
rotation curves are shown as solid lines, and +cr 1D envelopes are represented by dashed lines. Bottom: Histogram of amplitudes 
Q (left), y rot (centre) and A (right) from 10000 Monte Carlo simulations of RV distributions with no rotation. The optimal values 
for the observed RV configuration are indicated by dashed vertical lines. Confidence levels of 94.6%, 93.6% and 95.6% are found 
respectively for the best-fit amplitude of each model. 



density minima coinciding with the rotation axis (i.e. PAo~45° 
and 225°) and density maxima 90° away from the rotation axis 
(i.e. PAo~135 and 315°). However, only one minimum is seen at 
a position angle of ~ 100- 1 20° in the K s - and //-b and azimuthal 
density profiles of R136 bv lCampbell et alj d2010t ). 

Given the very young age of R136, the cause of the rotation 
needs to be looked for in the details of the formation process 
of the cluster. This is a short and complicated phase in which 
various physical processes, with their respective time-scales, op- 
erate simultaneously. A logical starting point is to see whether 
giant molecular clouds (GMCs), the birth sites of YMCs, ro- 
tate. iRosolowskv et alj d2003l) found that GMCs in M33 have 
non-zero angular momentum. From simple arguments based on 
differential rotation in a galactic potential and self gravity we 
expect that the rotation of GMCs should be prograde w ith the 
orbit in the galaxy. However, IRosolowskv et al. (l2003h found 
that 40% of the GMCs have retrograde motions, which sup- 
ports a scenario in which GMCs form through both agglom- 
eration and self gravity and the angul ar momentum i s the re- 
sult of the dumpiness of the gas dDobbs et al.L 1201 lb . As we 
pointed out above, a merger is anothe r way to give rise to rota- 
tion. Interestingly, Sabb i et ail d2012l) found a dual structure in 
the density of low mass stars in R136 that possibly hints at a 
relatively recent merger event of the main core of R136 and a 
second clump or cluster, but only about three of our targets are 
located in this second clump. 

To establish whether the rotation in old GCs is a remnant of 
their formation, we need to know if the angular momentum can 
survive for a Hubble time of dynamical evolutio n. During the 
evolution, angular momen tum is diffused outward ( Fa ll & Frenkl 
Il985t lEinsel & SpurzemL 1 19991) and ultimately lost through the 
escape of stars with high a ngular momentum dAgekianl Il958t 
Shap iro & MarchanH 119761) . This process operates on a relax- 
ation time and the angular momentum reduces after a fixed num- 



ber of elapsed relaxation times. The relaxation time of expanding 
clusters grows roughly linearly in time, which makes the number 
of elapsed relaxation times grow slowly, namely as a logarithm 
of the age. Because the majority of GCs are in this expansion 
phase dGieles et all 1201 lb . we may expect rotational signatures 
to still be present after a Hubble time because they have not 
evolved enough. For clusters that have entered the 'mass-loss' 
phase, we do not expect the rotation to survive. 

5. Conclusion 

We presented evidence that the young massive cluster R136 is 
rotating with a rotational velocity amplitude of about 3kms _1 , 
which implies that at least ~20% of its total kinetic energy is 
in rotation. Obviously, RV measurements of more stars in this 
cluster would be desirable to better populate the rotation curve 
and confirm the rotational signal with a confidence level higher 
than the current 95%. Given the young age of R136, our results 
suggests that star clusters may form with a significant amount of 
angular momentum. This will place useful constraints on mod- 
els of cluster formation. We finally argued that the rotation of 
globular clusters could originate from their formation, but this is 
clearly a topic where more detailed numerical investigations are 
welcome. 
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